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Instabilities for boson-fermion mixed condensates of trapped Alkali atoms due to the boson- 
fermion attractive interaction are studied using a variational method. Three regions are shown 
for their instabilities according to the boson-fermion interaction strength: stable, meta-stable and 
unstable ones. The stability condition is obtained analytically from the asymptotic expansion of 
the variational total energy. The life-time of metastable states is discussed for tunneling decay, and 
is estimated to be very long. It suggests that, except near the unstable border, meta-stable mixed 
condensate should be almost-stable against clusterizations. The critical border between meta-stable 
and unstable phases is calculated numerically and is shown to be consistent with the M0lmer scaling 
condition. 



Discoveries of the Bose-Einstein condensates (BEC) for Alkali atoms |l|J^ and the Fermi degeneracy for trapped 
atoms ^ encourage the study for boson-fermion mixed condensates of trapped atoms. Originally, the boson- 
fermion condensates have been studied in bulk systems such as ^He-^He and hydrogen-deuteron systems Studies 
of trapped mixed condensates have progressed recently; the static properties of 39.4i-|^_40j^ condensate (with repulsive 
boson-boson interaction) , dynamical expansions after the removal of trapping potentials Q , and the instabihty 
changes of ^Li-^Li system with attractive boson-boson interaction [ p^ . 

In this paper, we study the instabilities and collapses of the polarized boson-fermion mixed condensates due to 
the attractive boson-fermion interaction. In BEC, such instabilities are observed in several systems: 1) the trapped 
meta-stable ^Li BEC (observed experimentally |11 1) equilibrated between the attractive boson-boson interaction and 
the kinetic pressure due to the finite confinement |12[, 2) two-component uniform BEC (from bosons I and 2), whose 
stability condition is given by (711(722 > 5i2 (ffy • the coupling constant between bosons i and j). In the latter case, 
the stability condition depends only on the interaction strength, but the condition for the case 1) has also the particle 
number dependence fl^ . 

Let us consider the case of the trapped boson-fermion mixed condensate. In the present paper, we assume T = 
and a spherical harmonic-oscillator for the trapping potential; extension to the deformed potential is straightforward. 
Using the Thomas- Fermi approximation for the fermion degree of freedom (appropriate for large Nf condensate), 
the total energy of the system becomes a functional of the boson order-parameter $(r) and the fermion density 
distribution n/(r): 
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where oj and m are the angular frequency for the boson/fermion trapping potential and the boson/fermion mass 
(wf, — ujf = OJ and and nih = ruf = m are assumed for simplicity). In low-density gas, the boson-boson/boson- 
fermion coupling constants g, h are represented by the s-wave scattering lengths a;,f, and a{,/: g = Anh^ abb / rUh and 
h = 2i:h^ aiif /rrir where rrir = rrn,mf /{rrn, -t- to/). In eq. (Q), the fermion-fermion interaction has been neglected; the 
elastic fermion-fermion s-wave scattering for a polarized gas is absent because of the Pauli blocking effects and the 
p-wave scattering is suppressed below 100 /iK [Q. 

Evaluating the total energy with the variational method, we take the Gaussian ansatz for the boson order-parameter: 
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where x = r/^ is a radial distance scaled by the harmonic-oscillator length ^ = (Ti/tocj)^/^^ ^nd Nb = J (Px\^{x; R)\'^ 
is total boson number. The variational parameter R in corresponds to the root-mean-square (rms) radius of the 
boson density distribution. This kind of variational functions have also been used in the instability investigations of 



'^Li-BEC [|2j and ^Li-^Li system In the Thomas-Fermi approximation with (^. the fermion density distribution 
n/(r) becomes 

1 r 1 - 

nf{x; R) = -^ [nf{R) - - ^)|2j ^ (3) 

where h = 2h/ (huj^^). This fermion density vanishes at the turning point A which is determined by 

F{A) ^fif-K^- ^e-l^ = 0, (4) 

where C = {'i/2Tif/'^hNi,. We take n/(x; i?) ee when a; > A. 

The scaled chemical potential = 2^f/(h(jj) in eq. is determined by the normalization condition J (Pxnf{x) = 
Nf. As a result, the two parameters A and p,f become functions of R. 

Using eqs. ^ and (||), the total energy E{R) = E[^{x; R),nf{x; R)] becomes 



E{R)/{huo) = Nb 



+ Nf {he{R) + Iho{R) + CR-^hf{R)} , (5) 
where g = 2g/(huj^^) and integrals Ike,ho,bf{R) are defined by 



/fee(i?) = ^ r dxX^F{xf/\ 



/^o(i?) = :^ r dxx^F{xf/\ 
SnNf Jq 



hf{R) = r dxx^e-^F{xf/\ (6) 

One of the most important candidates for mixed condensate is the 39:4i]^_40j^ (boson-fermion) system. Their 
scattering lengths are not well fixed at present, and different values have been reported experimentally p5|-p^. In 
the present calculations, to estimate the qualitative stability behaviors, we take the value for af,f,(''^K) in ref. [T5| , 
and it gives g = 0.2 for the boson-boson interactions with w = 450 Hz. For the attractive boson-fermion interaction 
h, several negative values are taken: a = \h\/g = (0 ^ 3). It should be noted that the interaction strength can be 
shifted using the Feshbach resonance phenomena Q , whose applications for K-atom have been discussed in |^,^ . 

In fig. 1, numerical results for E{R) in (|) are plotted for a = \h\/g = 0.1, 1.0, 2.0, 2.5, 3 (lines a-e) with N = 10^ 
where three kinds of patterns can be read off: stable (a), meta-stable (b,c) and unstable states (d,e). In weak boson- 
fermion interaction, the system is stable against the variation of R, and has an absolute minimum as an equilibrium 
state. For the intermediate strength, it becomes meta-stable with one local minimum (lines b,c), and, in strong 
attractive interactions, the minimum disappears and the system become unstable (lines d,e). It should be noted that 
similar meta-stable states appear in the BEC with attractive interactions (^Li) ||l2|, but no stable states can exist in 
that case. 

As can be found in fig. 1 , the stability of mixed condensate can be judged from the small R behavior of E(R): positive 
divergence in small R suggests stable state. To obtain the stability condition, we consider analytically the asymptotic 
expansion of E{R) in i? ^ 1. For that purpose, the leading order terms of fif{R) and A(i?) in small R should 
be determined from eq. (Q) and the normalization condition: J d^xnf{x;R) — Nf. To evaluate them, we assume 
A{R)/R<t: 1 wheni?< 1, and expand the Gaussian function to the order of (x/i?)^: exp{-3x'^ /2R'^) ~ l-{3x'^ /2R'^). 
The consistency of these assumptions will be shown later, and can be also checked numerically. These assumptions 
makes the integral in the normalization condition evaluated analytically, and we obtain 

It should be noted that eq. (||) leads to A/R cx i?^/"* ^ 1 for i? <C 1, which shows the consistency of our assumption. 
Using eqs. (R) and (@), in the leading order in small R, the E{R) becomes 
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so that, in small R region, E{R) is dominated by i?^'^-term, and its positivity gives the stability condition: 

a^M<2-^/'^. (10) 
g Nf 

For the cases of fig. 1 , it becomes a < 0.18, which shows that only line (a) is for stable state. This condition is also 
plotted as the dashed line in fig. 2 . 

From the above stability condition, lines (b) and (c) in fig. 1 are found to be meta-stable, so that, in principle, the 
equilibrium states (at local minima R = Req) should collapse by the quantum tunneling effects into those of i? = 
through the potential barrier between R = and Req. To study the collapses of meta-stable states, we estimate 
their collective tunneling life-time applying the Gamow-theory of the nuclear a-decay to these states. For the 
order estimation of the tunneling life-time, as an approximation for E{R), we use the simplified linear-plus- harmonic 
oscillator potential V{R): 

..(J,. fVi{R) = FR-G, iO< R<Rt) .... 

^^^>-\V2(R)^^Am^{R-Re,f + Vo. {Rt<R) ^^^> 

As shown in appendix A, the tunneling life time for V{R) becomes t^^^ — £'exp(— M^), where Z? is a staying 
probability on barrier surface per unit time and exp(— M^) is a transmission coefficient for the potential barrier. The 
explicit formulae for D and exp(— M^) are 
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where Re is a WKB turning point for the energy E and Ay = Vi{Rt) — Vq is & barrier height (see fig. 3). ano 



■y/h/MU is har monic oscillator length around equilibrium point. The derivation of eqs. ( p^ , p^ is given in appendix 
A. 

To apply the above formula for the metastable condensates, we take M = Nm and = 2u). These values have been 
determined to be the values in harmonic oscillator potential around the equilibrium point in the non-interacting case. 
For parameters (i?OiK))j (-Rt,Vi(i?t)) and Re, we use the values that can be obtained from the numerically calculated 
E{R): the equilibrium, maximum, WKB turning points (for E = Vq + Ml/ 2). 

For the present case of = 10® and a = 2.0 we obtain Rt~ Re ^ 1, Ra — Rt ^ 3 in the unit of h/muj, and AV ~ 40 
in the unit of Nhu; (see the line (c) in fig. 1). Using these values for eqs. ( p^HT^ ), the staying probability and the 
transmission coefficient become D ^ lO^uj and W ^ 10^. As a result, we obtain r^t^/uj ~ 10^ x exp(— 10^), and the 
life-time Td becomes very longer than, for example, that of clusterization by many-body collisions, ~ (1 ~ 10) sec. 
Thus, meta-stable states should be really almost stable except the ones extremely close to the unstable border (such 
as line d). This long lifetime has originated from coUectiveness of very many particles, and is consistent with that 
for meta-stable EEC To shorten the tunneling life-time and make meta-stable states collapse, the parametric 

excitation with to ~ uj{t) may give an interesting possibility, which will be discussed in another paper p^ . 

For unstable states (e. g. line e in fig. 1), using eqs. (0,^), the leading i?-dependence of the fermion density at 
a; = and its rms radius are calculated: n/(0) (x R~~ and {xf)rms oc R^ , which are also checked numerically. They 
show that, in the collapses of metastable and unstable states, both boson and fermion distributions become localized 
and compressed, just like gravitational collapses of massive stars. Those collapses can also be found in the instability 
by the collective excitations evaluated with the sum-rule method p^ , where the boson- fermion in-phase monopole 
excitation becomes the zero mode. 

Finally, we comment about the critical condition for the unstable regions, which has originally been obtained by 
M0lmer Q with the Thomas-Fermi approximations for both fermion and boson distribution in form of the scaling 
relation: h'^/g cx N^p/'-^p+^^ for V cx r^-type trapping potential (the derivation of the scaling law by M0lmer, see 
appendix B). For the harmonic oscillator potential {p — 2), it becomes a'^N^^ ~ const/y. In the present framework 
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of the variational method, the critical a between the meta-stable and unstable states have been evaluated for several 
N and plotted in fig. 2 (filled circles). From them, we can find that the variational results are also consistent with 
the M0lmer scaling law: the solid line added in fig. 2 in the unstable region. 

In summary, for stability of boson-fermion mixed condensates with boson-boson repulsive and boson-fcrmion attrac- 
tive interactions, the present variational calculations give the phase diagram with three regions: stable, meta-stable 
and unstable ones. They are clearly shown in fig. 2 in the case of Nf — iVf, = N and ujf = uji, = lu. To estimate 
borders of stability regions generally, we studied their critical conditions analytically, from which we can find that 
the diagram should have qualitatively the similar pattern also in general cases of the mixed condensates. Finally, it 
should be noted that, in the boson-fermion collapsing processes of unstable and meta-stable states, the increasing 
boson/fermion densities make the two-body short-range interactions (e.g. p-wave scattering processes) and three-body 
(or higher) collisions more effective. These interactions, which are not essential for low-density condensates, may play 
important roles through the fermion-paired superfluid formations and clusterization into metal states by tunneling 
phenomena ||2^ . They will be discussed in a future publication . 



APPENDIX A: LIFE TIME OF COLLECTIVE TUNNELING EFFECT IN METASTABLE REGION 



To estimate a life time of metastable mixed condensates, we consider the collective tunneling effect with the collective 
variable R (the boson radius) and its effective energy E{R) as a collective potential. The effective energy E{R) is 
obtained from (^, and, for the metastable condensates, a typical shape of it is the line (c) in Fig. 1. We approximate 
this E{R) with the linear-plus-harmonic-oscillator type potential eq. ( pl| ) (Fig. 3). M is inertia mass for the collective 
variable R. The meanings of other parameters {F, G, fl, Req, Rt) can be read off in Fig. 3, and they are fixed in order 
that the potential V{R) reproduce the E{R) numerically obtained by (||). This approximation should be enough for 
the order estimation of the life time. For the kinetic term, we take 



2M dR'^ ' 



(Al) 



In V{R)^ the metastable state ipMix) before tunneling is approximately given by the ground-state wave function in 
the harmonic oscillator potential V2 {R) : 



^Af(i?) 



: exp 



{R — Re 



''HO 



(A2) 



where qho = ^/h/MJi is a harmonic oscillator length. This state has the zero-point energy measured from Vq. 

The continuum state ipniR) after the tunneling decay is obtained by the wave function in the linear potential Vi(i?), 
and its Schrodinger equation becomes 



2M di?2 



E ~ {FR ~ G) 



^d{R) = 0. 



(A3) 



To solve the above equation, we use the WKB approximation. It should be noticed that, because the state energy 
is lower than the potential maximum {E < FRt — G), a turning point exist at Re — {E + G)/F. Thus, the WKB 
connection formula should be used at i?^ for the continuum state ipo- 
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where A = ^2M/(TTh^) and nK{R) = ^2MF\Re - R\. 

Let us consider the tunneling decay rate F from iPm{R) to ^je){R). Using the golden rule of Fermi, it becomes 



F = 



27r 



dRi:*jj{R)[V{R) - V2{RMm{R) 
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where Em = Vo + hn/2 is the energy of Because V(i?) - V2{R) = tor Rt < R and Vm(0), -0^/(0) < 1, eq. 

( [A^ ) becomes 
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A4, A5), the tunneUng rate T is shown to have the form F = Z?exp(— M^) : £> can be interpreted as 



Using eqs. (| 

the staying probabihty Per unit time and exp(— M^) is the transition coefficient. The exphcit formulae for D and W 
are 
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where AV = Vi{Rt) — Vb is a barrier height. To represent the parameters F and G in Vi by AT^, we have used 



Vi{Rt) = FRt - G, Vi{Re) = FRe -G = E = Vo + ^fiil. 



(AlO) 



The tunnehng hfe-time is defined as an inverse of decay rate: Tc; 



APPENDIX B: CRITICAL CONDITION FOR UNSTABLE REGION 



In this appendix, we rederive the M0lmer scahng law which gives the critical condition for the unstable region 
[^. We assume the Thomas- Fermi approximations both for the boson and fermion density distributions. In that 
approximations, the density distributions are given as solutions of the Thomas-Fermi equations: 



gnb{x) + + hnf{x) ~ fib, 
[67r^n/(a;)]t + x'^ + hnb{x) = 



(Bl) 
(B2) 



where nb{x) is the boson density distribution scaled by the harmonic-oscillator length ^ — {H/mujY^'^ inb{x) = 
|$(a;)p^^). The scaled chemical potentials, fibj, and coupling constants, g and h, have been defined in the main body 
of this paper. 



Eliminating nb{x) in ( p32D by (Bl), we obtain the equation F[nf{x)] — G[nf{x)] to determine the fermion density 
nf{x), where 



\h\ 



F[nf{x)] = [67r^n/(x)] 3 , G[nf{x)] = fif + y Mb - (^1 + y j ^ + ^^fi^)- 
We concentrate on the central density ?t./(0). In order to determine nf{0) two conditions should be satisfied: 

FlufiO)] = GiufiO)], |^[n/(0)] > |^[n/(0)], 



Evaluating equations in (B4) with (B3), we obtain the critical condition for the unstable region: 



(B3) 



(B4) 



(B5) 



To evaluate the right-hand side of (p^), we should use the relations between jlbj and the boson/fermion particle 
number Nbj, which are obtained by solving the Thomas- Fermi equations in (B2). Here we assume that jib = and 
flf = 2{6Nfy/^ {flf for a free fermion system) for (B£). Consequently, we obtain the M0lmer scaling relation: 



g 61/6^6 ^ ^ 



(B6) 
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It should be noted that the coefRcient 12.0 in (B6) are close to the value 13.8 which is obtained by M0lmer j|] with 
numerically evaluating the Thomas-Fermi equations. 
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FIG. 1. The total energy variation E[R) / (Nhw) of the boson-fermion mixed condensate with the Gaussian boson distribution 
$(1; R) of radius R: E{R) and $(a;; R) are defined in The Thomas-Fermi density function (^) is applied for the fermion 

distribution. N = Nt = Nf = 10^, g = 0.2, and a = h/g = 0.1, 1.0, 2, 2.5, 3 for lines a-e. 

FIG. 2. Stability phase diagram of the boson-fermion mixed condensate for N = Ni, = Nf and a = h/g when g = 0.2. Three 
regions exist in it: stable (S), meta-stable (MS) and unstable (US) ones. The dashed line at a = 0.18 corresponds to the border 
of stable region, and the solid line shows the M0lmer scaling law between meta-stable and unstable regions. The open circles 
correspond to the parameters of the states in fig. 1, and the filled circles show numerically confirmed critical states between 
meta-stable and unstable regions. 

FIG. 3. Simplified potential for meta-stable condensates. The Vi(-R) and V2{R) are the linear and harmonic oscillator parts 
of the potential V^(-R). The Req, Rt, Re are the equilibrium, maximum and the WKB turning points. Vb is the equilibrium 
energy of the potential. 
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